Remote Sensing of Environment 114 (2010) 1805-1816 



ELSEVIER 


Contents lists available at ScienceDirect 

Remote Sensing of Environment 

journal homepage: www.elsevier.com/locate/rse 



Land surface phenology from MODIS: Characterization of the Collection 5 global land 
cover dynamics product 

Sangram Ganguly a ’*, Mark A. Friedl b , Bin Tan c , Xiaoyang Zhang d , Manish Verma b 

a Bay Area Environmental Research Institute/NASA Ames Research Center, MS 242-4, Moffett Field, CA 94035, USA 
h Department of Geography & Environment, Boston University, 675 Commonwealth Avenue, Boston, MA 02215, USA 
c Earth Resources Technology, Inc., NASA Goddard Space Flight Center, Code 614.5, Greenbelt, MD 20771, USA 
d Earth Resources Technology, Inc., 10810 Guilford Road, Suite 105, Annapolis Junction, MD 20701, USA 


ARTICLE INFO 


ABSTRACT 


Article history: 

Received 29 September 2009 
Received in revised form 8 April 2010 
Accepted 10 April 2010 


Keywords: 

MODIS phenology 
Land cover dynamics 
Land cover 
Vegetation dynamics 
MODIS Collection 5 


Information related to land surface phenology is important for a variety of applications. For example, 
phenology is widely used as a diagnostic of ecosystem response to global change. In addition, phenology 
influences seasonal scale fluxes of water, energy, and carbon between the land surface and atmosphere. 
Increasingly, the importance of phenology for studies of habitat and biodiversity is also being recognized. 
While many data sets related to plant phenology have been collected at specific sites or in networks focused 
on individual plants or plant species, remote sensing provides the only way to observe and monitor 
phenology over large scales and at regular intervals. The MODIS Global Land Cover Dynamics Product was 
developed to support investigations that require regional to global scale information related to spatio- 
temporal dynamics in land surface phenology. Here we describe the Collection 5 version of this product, 
which represents a substantial refinement relative to the Collection 4 product. This new version provides 
information related to land surface phenology at higher spatial resolution than Collection 4 (500-m vs. 1-km), 
and is based on 8-day instead of 16-day input data. The paper presents a brief overview of the algorithm, 
followed by an assessment of the product. To this end, we present (1) a comparison of results from 
Collection 5 versus Collection 4 for selected MODIS tiles that span a range of climate and ecological 
conditions, (2) a characterization of interannual variation in Collections 4 and 5 data for North America 
from 2001 to 2006, and (3) a comparison of Collection 5 results against ground observations for two forest 
sites in the northeastern United States. Results show that the Collection 5 product is qualitatively similar to 
Collection 4. However, Collection 5 has fewer missing values outside of regions with persistent cloud cover 
and atmospheric aerosols. Interannual variability in Collection 5 is consistent with expected ranges of 
variance suggesting that the algorithm is reliable and robust, except in the tropics where some systematic 
differences are observed. Finally, comparisons with ground data suggest that the algorithm is performing 
well, but that end of season metrics associated with vegetation senescence and dormancy have higher 
uncertainties than start of season metrics. 

© 2010 Elsevier Inc. All rights reserved. 


1. Introduction 

Investigations focused on monitoring and modeling biospheric 
processes require accurate information related to spatio-temporal 
dynamics in ecosystem properties. Because vegetation phenology 
affects terrestrial carbon cycling across a wide range of ecosystem and 
climate regimes (Baldocchi et al., 2001; Churkina et al., 2005; 
Richarson et al., 2009), accurate information related to phenology is 
important to studies of regional-to-global carbon budgets. The 
presence of leaves also influences land surface albedo (Moore et al., 
1996; Ollinger et al., 2008) and exerts strong control on surface 
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radiation budgets and the partitioning of net radiation between latent 
and sensible heat fluxes (Chen & Dudhia, 2001; Yang et al., 2001). 
Thus, the phenological dynamics of vegetated ecosystems influence a 
host of eco-physiological processes that affect hydrologic processes 
(Hogg et al., 2000), nutrient-cycling, (Cooke & Weih, 2005), and land- 
atmosphere interactions (Heimann et al., 1998). 

In recent years, growing season dynamics including shifts in the 
timing of bud burst, leaf development, senescence, and changes in 
growing season length have been widely studied in the context of 
ecosystem responses to climate change (Cleland et al., 2007). Studies 
using AVHRR data have concluded that northern hemisphere 
temperate and boreal regions (~40°-70° N) experienced increased 
growing season greenness related to surface warming during the 
period 1981 to 1999 (Myneni et al., 1997; Zhou et al., 2001). More 
recent studies utilizing a longer record of AVHRR data suggest a more 
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complex pattern with evidence of “browning” trends in the boreal 
forests of Southern Alaska, Canada, and in the interior of Russia 
(Angert et al, 2005; Ganguly et al., 2008; Goetz et al., 2005; Zhang et 
al., 2007). Complex phenological responses have also been observed 
in controlled experiments where warming accelerated greening of 
plant canopies, but elevated C0 2 and nitrogen fertilization delayed 
flowering (Cleland et al., 2007). These biophysical and biochemical 
processes both influence and are diagnostic of ecosystem-climate 
interactions. As a consequence, there is substantial need to accurately 
characterize the phenology of ecosystems, and by extension, the 
response of ecosystems to changes in climate (Morisette et al., 2009). 

Moderate resolution satellite remote sensing provides global high 
temporal frequency measurements of land surface properties and is 
therefore well suited for monitoring seasonal-to-decadal patterns and 
trends in regional-to-global phenology (de Beurs & Henebry, 2005; 
Reed et al., 1994; White et al., 1997; Zhang et al., 2003). Landsat 
MSS was the first space-borne sensor used to characterize the sea- 
sonality of vegetation at landscape and regional scales (Thompson & 
Wehmanen, 1979). However, detecting phenological transition dates 
requires higher temporal resolution than is afforded by Landsat-class 
instruments, and coarse-to-moderate spatial resolution sensors such 
as the AVHRR (Goward et al., 1985), MODIS (Zhang et al., 2003) and 
SPOT-VEGETATION (Delbart et al., 2006) are more commonly used for 
this purpose. Indeed, the utility of such sensors for studies of land 
surface phenology has been established for over 20 years (Justice et 
al., 1 985). Over the last two decades, a number of different approaches 
have been developed for this purpose. The most common of these 
methods include threshold-based techniques (Jonsson & Eklundh, 
2002; White et al., 1997), methods based on spectral analysis 
(Jakubauskas et al., 2001; Moody & Johnson, 2001), and inflection 
point estimation in time series of vegetation indices (Moulin et al., 
1997; Zhang et al., 2003). All methods use time series of vegetation 
indices to identify the timing of phenological transition dates such as 
the start and end of the growing season. 

Since 2000, the Moderate Resolution Imaging Spectroradiometer 
(MODIS) has provided an excellent basis for regional-to-global scale 
studies of land surface phenology (Ahl et al., 2006; Fisher & Mustard, 
2007; Zhang et al., 2003, 2006). The objective of this paper is to 
present an overview and characterization of the new Collection 5 (C5) 
MODIS Global Land Cover Dynamics (MLCD) product, which is pro- 
duced globally at a spatial resolution of 500-m and is available from 
2001 to present. Below we describe refinements to the MLCD algo- 
rithm and present a characterization of the C5 product via ( 1 ) com- 
parison with Collection 4 (C4) data, (2) assessment of variability in C5 
results across multiple years, and (3) comparison of C4 and C5 data 
with field data. 


2. MODIS land cover dynamics algorithm and product 

2.1. Product overview , input data, and pre-processing 

The MLCD algorithm, as presented in Zhang et al. (2003, 2006), is 
used for this work. This algorithm characterizes vegetation growth 
cycles using four transition dates estimated from time series of MODIS 
enhanced vegetation index (EVI) data: (1) greenup: the date of on- 
set of EVI increase; (2) maturity: the date of onset of EVI maximum; 
(3) senescence: the date of onset of EVI decrease; and (4) dormancy: 
the date of onset of EVI minimum. These transition dates correspond 
to the first four science data sets (SDS) associated with the MLCD 
product (Table 1). Note that the SDS names use “greenness”, not EVI. 
To be more precise, and because these transition dates are estimated 
from temporal dynamics in the EVI, here we specifically refer to 
transitions in EVI (onset of increase in EVI, etc.). In Section 3.4 we 
compare these dates with measurements of forest canopy phenology 
collected on the ground. 


Table 1 

Science data sets included in the MODIS Land Cover Dynamics Product (combined Terra 
and Aqua Product Identifier MCD12Q2). 


Science data set name 

Units 

Data 

format 

Valid range 

Scale 

factor 

Onset_Greenness_Increase 

Onset_Greenness_Maximum 

Onset_Greenness_Decrease 

Onset_Greenness_Minimum 

Day of 
year 

16 bit 
signed 
integer 

[1366] 

N/A 

EVI_Onset_Greenness_Increase 

EVI_Onset_Greenness_Maximum 

EVI units 


[-10,000,10,000] 

10 4 

EVI_Growing_Season_Area 

XEVI 


[-3660,3660] 

10 2 

Phenology_Quality 

N/A 

8 bit 
signed 

Bit fields 

N/A 


To illustrate a potential application of these data, Fig. 1 shows the 
mean growing season length for North America for the period 2001- 
2006 derived from MLCD data. To generate this figure, the growing 
season length for each pixel in each year was calculated as the 
difference between the time of onset of EVI increase and the time of 
onset of EVI minimum. Because growing season length is an important 
variable that affects the time period during which photosynthesis 
takes place, information derived from the MLCD product is useful for 
understanding how seasonal and interannual climate forcing affects 
plant growth, and by extension, gross and net ecosystem exchange 
over large scales. Recently, Medvigy et al. (2009) and Richardson et al. 
(in press) have demonstrated the utility of MLCD data for model- 
based and empirical studies of seasonal-scale carbon budgets, 
respectively. In Section 3.3 we present sample results quantifying 
the nature and magnitude of interannual variability in growing season 
length for North America based on the MLCD product. 

The C5 MLCD product is distributed as separate science data sets 
(SDS) within a single file stored in Hierarchical Data Format. All SDS 
layers are produced at 500-m spatial resolution in the MODIS 
sinusoidal projection. The product includes eight SDSs, which are 
summarized in Table 1. Note that vegetation can have more than one 
growth cycle during any 12-month period. Each SDS field therefore 
includes two 16-bit values, one for each of two possible cycles in a 12- 
month period. To accommodate differences in northern and southern 
hemisphere seasonality, the MLCD product is created at six-month 
time steps. 

As we indicated above, vegetation growth cycles are estimated 
using MODIS EVI data. However, rather than using the standard MODIS 
EVI product (Huete et al., 2002), the MLCD product is generated using 
EVI computed from MODIS nadir bidirectional reflectance distribu- 
tion function (BRDF)-adjusted reflectance (NBAR) data (Schaaf et al., 
2002). In Collection 5, NBAR data is produced at 500-m spatial resolu- 
tion and 8-day time steps based on overlapping 16-day periods. The 
effects of clouds, variable view angles, and atmospheric aerosols 
are minimized in NBAR data. The EVI is used because it provides 
greater dynamic range than the normalized difference vegetation 
index (Huete et al., 2002). 

The NBAR product also provides a snow and ice flag in its quality 
assurance field, indicating whether the data are acquired over a snow- 
covered or snow-free surface. This information is critical to our ap- 
proach, as the presence of snow lowers EVI values during the winter 
months (Delbart et al., 2006; Dye & Tucker, 2003). In addition, we 
utilize the 1-km C5 MODIS land surface temperature (LST) product, 
also available at 8-day time intervals (Wan et al., 2002). At each pixel 
and time step, EVI values identified to contain snow are replaced with 
a background EVI value, which is defined as the most recent snow-free 
EVI value. Gaps caused by clouds and atmospheric aerosols are filled 
using a centered three-date moving window average. As a final step, 
the EVI time series are smoothed using a local median-value moving- 
window technique (Zhang et al., 2006). This last step has the effect of 
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smoothing the time series by removing EVI values that depart sub- 
stantially from the local trend in the time series. 

2.2. Algorithm overview 

The MLCD algorithm is described in Zhang et al. (2006). Here we 
provide a summary of the algorithm. Periods with sustained EVI 
increase and decrease are determined using moving windows 
composed of five 8-day NBAR EVI values. Because small decreases or 
increases can be caused by local or transient processes unrelated 
to growth cycles, two heuristics are used to exclude such variation: 

( 1 ) the change in EVI within any period of EVI increase or decrease 
must be larger than 35% of the annual range in EVI for that pixel; and 

(2) the ratio of the local maximum EVI to the annual maximum EVI 
should be at least 0.7. This approach screens out short-term variation 
in EVI unrelated to growth or senescence cycles (Zhang et al., 2006). 

Following Zhang et al. (2003), temporal variation in EVI for a single 
growth or senescence cycle is modeled using a sigmoid function of the 
form: 

m(t)= 1+e ' + ^ + d (D 

Where t is time in days, a and b are fitted parameters that are 
estimated using a non-linear least squares algorithm, c is the am- 
plitude of EVI variation, and d is the minimum snow-free EVI value at 
each pixel. The latter value is defined as the minimum stable value 
without snow present in the pixel. Using this framework, phenolog- 
ical transition dates are identified based on the rate of change in 
curvature of the fitted logistic models at each pixel (Zhang et al., 2003, 
2006). In addition, the algorithm computes the growing season 
integrated NBAR EVI (the sum of the modeled daily EVI values from 
the onset of EVI increase to the onset of EVI minimum), and records 
the EVI values corresponding to the transition dates in the modeled 
EVI time series. Key differences relative to the C4 product are that ( 1 ) 
C5 is produced at 500-m spatial resolution based on 8-day instead of 
16-day periods; (2) the input EVI data benefit from refinements to 
upstream products (surface reflectance and NBAR data); and (3) the 
algorithm no longer applies a 3 by 3 moving window average at each 
pixel to remove high frequency spatial variance. Thus, the spatial 


resolution of the product has increased from an effective spatial 
resolution of 3-km to 500-m. 

3. Analysis and results 

In this section we present results from four different sets of 
analysis. Because a global comparison is not practical, we present 
results for six MODIS tiles (h23v02, hl2v03, hl2v04, h27v05, h20v07, 
and h20v08) encompassing a wide range of climate and land cover 
types (Fig. 2). Hereafter, we refer to these as the subarctic Eurasian, 
boreal North American, temperate North American, temperate Asian, 
arid tropical, and humid tropical tiles, respectively. These tiles include 
tundra, boreal, temperate, sub-tropical and tropical climate regimes, 
and encompass biome types that span the entire range of temperature 
and moisture-limited life forms. Two of these tiles include regions 
with substantial human management of land cover and land use. All of 
the IGBP land cover types are included in this sample. With the 
exception of the IGBP closed shrublands (which is a small class), 
water, urban, and permanent snow and ice classes, each IGBP class is 
included as one of three dominant classes in one or more of these tiles 
(Table 2). 

We use these tiles throughout the analysis presented in Sections 3.1, 
3.2, and 3.3. We begin in Section 3.1 by characterizing the nature and 
magnitude of missing data in the C5 MLCD product. In Section 3.2 we 
present a comparison of the MLCD product for C4 versus C5. The goal of 
this analysis is to provide a qualitative and quantitative assessment of 
differences and similarities between the two versions of the product. In 
Section 3.3, we examine C5 MLCD results across different years and 
assess the nature and magnitude of interannual variation in C5 data sets. 
Finally, in Section 3.4 we compare results from the C4 and C5 products 
against ground data from two sites, which provides a quantitative 
assessment of the MLCD product quality relative to phenology 
measured on the ground. In Sections 3.2 and 3.3 the 1-km C4 land 
cover dynamics product is used as a benchmark for comparison with the 
C5 product. To allow this comparison, C4 data were resampled to 500-m 
using bi-cubic interpolation. The phenological metrics produced by the 
MLCD algorithm are the same for both C4 and C5, and aside from 
different input data and spatial resolutions, are directly comparable. The 
MODIS C5 global land cover type map (MCD12Q1 ; Friedl et al., 2010) is 
also used in this exercise to stratify the data into different land cover 
classes. 
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Fig. 2. Global sinusoidal MODIS tile map displaying the six tiles (solid black boxes) that are used in the analysis. The labels classify these tiles according to the dominant land cover 
type in each tile from the MODIS Land Cover Type product. 


3.1. Missing data 

We first explore the nature and magnitude of missing input data in 
C5 MLCD results relative to C4. The C4 product was produced for 2001 
through 2004 and the C5 product spans the time period from 2001 to 
present. Here we consider data from 2002 from each collection. Before 
addressing the more general question of missing values arising from 
clouds, atmospheric aerosols, and low illumination conditions we first 
discuss a specific issue. In 2001, instrument problems on the Terra- 
MODIS sensor from day-of-year (DOY) 166 (June 15) to DOY 183 (July 
2) resulted in a data gap during a critical period of the northern 
hemisphere growing season. As a consequence, the quality of MLCD 
retrievals during this period was significantly compromised in C4. In 
C5, MLCD input data for DOY 166-183 in 2001 were imputed using 
mean NBAR data for the same period from 2002 to 2006 at each pixel. 
Fig. 3 shows the original EVI time series and the gap-filled data for a 


pixel in central Massachusetts and illustrates the effectiveness of this 
solution. Prior to applying this adjustment, noticeable artifacts were 
present in MLCD results throughout the northern hemisphere. Visual 
inspection of C5 results confirms that this strategy removes the 
majority of these artifacts, which were also present in C4 data. That 
said, MLCD values during this period are necessarily of modestly 
lower quality relative to other time periods. 

More generally, missing values caused by clouds and aerosols are a 
major source of uncertainty in MLCD results (Zhang et al., 2009). 
However, because the 8-day C5 MODIS reflectance data provide 
higher frequency data with better cloud screening and atmospheric 
correction relative to C4, the C5 product has fewer missing values 
relative to C4. To illustrate, Fig. 4 shows the proportion of missing 
retrievals for greenup onset dates in C5 and C4 in 2002 for each of the 
six tiles shown in Fig. 2. Fig. 4 shows that C5 data have substantially 
fewer missing retrievals for the subarctic Eurasian, boreal North 


Table 2 

Summary statistics (mean, standard deviation) for MLCD SDS’s in each of the selected tiles, stratified by the three dominant land cover types in each. LC types are IGBP land cover 
classes from the MCD12Q1 product: 1. Evergreen needleleaf forest; 2: evergreen broadleaf forest; 3 deciduous needleleaf forest; 4: deciduous broadleaf forest; 5: mixed forest; 
7: open shrublands; 8: woody savannas; 9 savannas; 10: grasslands; 11: wetlands; 12 agriculture; 14: agricultural mosaic; 16: barren and sparsely vegetated. 


Tile ID 

LC 

type 

Pixels 

(xlO 4 ) 

% LC 

Onset EVI 
Increase 

Onset EVI 
Maximum 

Onset EVI 
Decrease 

Onset EVI 
Minimum 

Integrated EVI 

Maximum EVI 

Minimum EVI 

X 

<J 

X 

( 7 

X 

< 7 

X 

< 7 

X 

( 7 

X 

( 7 

X 

o 

Subarctic Eurasia 

3 

241.01 

43 

141 

5 

179 

5 

219 

3 

255 

3 

38.86 

1.14 

0.42 

0.0038 

0.20 

0.0060 


7 

194.96 

35 

157 

5 

189 

4 

219 

4 

246 

2 

28.18 

0.76 

0.36 

0.0036 

0.23 

0.0057 


8 

89.90 

16 

151 

5 

187 

4 

220 

4 

249 

2 

28.09 

0.97 

0.34 

0.0046 

0.19 

0.0058 

Boreal North America 

1 

255.90 

48 

135 

11 

188 

7 

219 

6 

268 

3 

39.75 

3.51 

0.35 

0.0046 

0.22 

0.0062 


11 

100.82 

19 

136 

12 

192 

7 

224 

5 

271 

3 

42.26 

4.44 

0.37 

0.0053 

0.23 

0.0060 


4 

55.59 

10 

138 

12 

192 

8 

223 

6 

269 

2 

42.02 

4.15 

0.38 

0.0057 

0.23 

0.0048 

Temperate North America 

5 

168.07 

40 

130 

4 

177 

5 

221 

2 

286 

4 

68.75 

2.12 

0.55 

0.0059 

0.27 

0.0097 


4 

92.64 

22 

117 

3 

163 

4 

223 

1 

304 

3 

100.15 

3.28 

0.70 

0.0132 

0.27 

0.0118 


14 

79.62 

19 

109 

3 

171 

2 

223 

3 

303 

2 

93.09 

2.72 

0.62 

0.0114 

0.25 

0.0091 

Arid Tropical 

16 

299.15 

52 

212 

11 

230 

9 

249 

7 

254 

16 

9.67 

1.07 

0.14 

0.0100 

0.09 

0.0009 


10 

164.67 

29 

182 

9 

232 

6 

255 

6 

287 

9 

32.79 

2.62 

0.33 

0.0136 

0.15 

0.0041 


7 

55.34 

9 

186 

6 

230 

6 

252 

5 

277 

8 

30.69 

2.27 

0.30 

0.0091 

0.14 

0.0028 

Humid Tropical 

2 

225.98 

39 

95 

14 

138 

12 

178 

26 

183 

33 

89.66 

26.13 

0.58 

0.0180 

0.30 

0.0972 


8 

202.73 

35 

60 

5 

186 

3 

251 

6 

276 

42 

124.47 

4.04 

0.58 

0.0050 

0.26 

0.0101 


9 

126.95 

22 

106 

6 

211 

3 

252 

2 

307 

16 

92.64 

2.15 

0.57 

0.0048 

0.20 

0.0098 

Temperate Asia 

12 

267.03 

59 

120 

17 

163 

13 

188 

10 

229 

7 

38.37 

4.90 

0.46 

0.0162 

0.11 

0.0351 


5 

58.96 

13 

92 

11 

155 

2 

209 

3 

277 

9 

82.75 

6.05 

0.59 

0.0161 

0.17 

0.0183 


14 

47.72 

10 

113 

11 

168 

7 

207 

8 

250 

14 

55.38 

3.81 

0.52 

0.0207 

0.12 

0.0158 
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Fig. 3. Sample time series for a pixel in central Massachusetts showing the original EVI 
time series for 2001 including the data gap in late June and early July (blue), and the 
2002-2006 mean EVI for the same period used to fill the missing data. Note that only 
data from the missing period were replaced. The rest of the time series is shown for 
completeness. 


American, temperate North American, and temperate Asian tiles. For 
the other two tiles the proportion of missing data is comparable 
(-40-50%), but is slightly higher in C5 relative to C4. This result 
reflects the fact that areas in the tropics and Asia are characterized by 
persistent cloud cover, high levels of atmospheric aerosols, or weak 
seasonality, all of which present substantial challenges for land 
surface phenology algorithms. To avoid generating spurious results, 
the MLCD algorithm is designed to be conservative and does not 
produce a result if data are missing during transition periods or if the 
amplitude in seasonal EVI is very low. This results in high levels of 
missing values in parts of the tropics and subtopics, and in areas with 
persistent high aerosol loading. In the specific case of the tiles 
examined here, land cover in the humid tropical and arid tropical tiles 
includes 39% evergreen broadleaf forest and 52% barren and sparsely 
vegetated land cover, respectively. Thus, the relatively high percent- 
age of missing values in these tiles is neither surprising nor 
inappropriate. In the temperate Asian tile, there are fewer missing 
values in C5, but the overall rate is still high. Land cover in this tile is 



Tile Identifiers 

Fig. 4. Proportion of missing retrievals for onset of EVI increase dates from the C5 and C4 
products for 2002. 


dominated by agriculture and mixed forests and the high percentages 
of missing values is related to persistent cloud cover and atmospheric 
aerosols in this region of Asia. 

3.2. Comparison with Collection 4 MLCD data 

We now present a comparison of the MLCD products for C5 and C4. 
Fig. 5 presents images showing retrieved transition dates for C5 (year 
2002) corresponding to the onset time of EVI increase, maximum, 
decrease, and minimum for the six MODIS tiles used in our analysis. To 
save space, here we present results from 2002 only. To assess the 
representativeness of 2002 results, we performed parallel analyses for 
2003-2006. No substantial differences were observed across years 
during this time period and we therefore concluded that presentation 
of results for 2002 is sufficient for the goals of this comparison. 

Fig. 6 shows frequency distributions for differences in transition 
dates for 2002 (C5 minus C4) for the six MODIS tiles identified above. 
Only pixels with valid retrievals in both data sets were included. 
Differences between C4 and C5 are generally within ±50 days. 
Indeed, most differences were of much smaller magnitude, and some 
differences are introduced by resampling the 1-km C4 data to 500-m 
spatial resolution. However, the arid subtropical and humid tropical 
tiles show substantial differences in the timing of onset of EVI 
increase. Specifically, while the primary mode of the distribution is 
close to zero, there is a clear secondary mode in both tiles associated 
with later greenup dates in C5 compared to C4. Unfortunately, C4 
NBAR data are no longer available for comparison. Fig. 7 presents the 
mean C5 EVI profile (±1 standard deviation) for pixels in the humid 
tropical tile where the difference in timing is between 30 and 70 days. 
This figure shows that the C5 result is reasonable, and that the 
algorithm may have been spuriously detecting early season greening 
in C4. 

Fig. 8 presents scatterplots showing mean transition dates for IGBP 
land cover types in C5 versus C4 for 2002. The class means are in close 
agreement for most tiles. The humid tropical tile shows substantial 
disagreement in the average date of onset for EVI decrease and EVI 
minimum for several classes (evergreen broadleaf and deciduous 
broadleaf forests, woody savannas, and mixed forests). Overall, how- 
ever, the comparisons presented in Figs. 6 and 8 indicate that despite 
differences in the spatial and temporal resolution in each collection, 
results from each of the products are generally quite consistent, with 
some explainable discrepancies in the tropics. 

3.3. Interannual variability with land cover 

In this section, we characterize the nature and magnitude of 
variability in C5 phenological metrics at interannual time scales from 
2001 to 2006. Table 2 presents summary statistics where, for each of 
the three most common classes in each tile, we present the overall 
mean value for each SDS listed in Table 1 (excluding phenology 
quality), along with corresponding standard deviations. The patterns 
in average phenology revealed in this table are consistent with 
expected biogeographic patterns: the timing and length of the 
growing season varies with climate and land cover, with unmanaged 
temperate and humid tropical land cover types showing longer 
growing seasons and higher overall greenness than high-latitude, arid 
tropical, and managed land cover types. 

Closer inspection of Table 2 reveals a number of patterns that merit 
further attention. Specifically, depending on geographic location, land 
cover exerts strong control on phenology. For example, the onset of 
EVI increase for IGBP classes 3 (subarctic Eurasia), 5 (temperate North 
America), 16 (arid subtropical), and 8 (humid tropical) are quite 
different from neighboring land cover types in these tiles. The onset of 
EVI minimum also shows substantial dependence with land cover 
(e.g., class 16 in the arid subtropical tile and all three classes for the 
humid subtropical and temperate Asian tiles). With relatively few 
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Fig. 5. Images showing the timing of onset for EVI increase, maximum, decrease and minimum for the six tiles used in the analysis (year 2002). Fig. 2 shows the location for each of 
these tiles. 
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Fig. 6. Frequency distributions for differences between the C5 and C4 Land Cover Dynamics products at 500-m spatial resolution. Note that the distributions have been normalized to 
scale from 0 to 1. 
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Fig. 7. Mean (± 1 standard deviation) EVI time series for pixels in tile h20v08 where C5 
greenup onset is between 30 and 70 days after C4 greenup. 


exceptions, average timing in the onset of EVI maximum and the onset 
of EVI decrease are generally quite uniform within classes. 

Equally important, the magnitude of interannual variability (i.e., 
standard deviations) in transition dates varies strongly with land 





cover. Year-to-year variability is relatively high for the onset of EVI 
increase in the boreal North American and temperate Asian tiles, and 
for the onset of EVI minimum in the humid tropical tile. Individual 
classes in specific tiles (e.g., class 2 in the humid tropical tile) also 
show large standard deviations. High variance in the onset of EVI 
increase in the temperate Asian tile is largely associated with 
agricultural land use. In the case of the North American boreal tile, 
which is dominated by the boreal forest biome, variability introduced 
by fires and regional “browning” trends may partly explain this 
pattern (Bunn & Goetz, 2006; Soja et al., 2007). More generally, large 
variance in the timing of phenology corresponds to classes and 
locations that present the greatest challenges for remote sensing: 
evergreen systems (classes 1, 2), locations with high percentages of 
missing input data (i.e., the temperate Asian and humid tropical tiles), 
and locations at high latitudes with low solar zenith angles and 
heterogeneous landscapes (boreal North America). Surprisingly, the 
subarctic Eurasian tile shows relatively low levels of variability in the 
timing of phenology within classes, which probably reflects the strong 
control that synoptic temperature patterns exert on phenology in this 
region. 

Variability in integrated EVI for class 2 in the humid tropical tile is 
an order of magnitude larger than for any other class or tile. This result 
reflects difficulties involved in remote sensing of phenology for 
evergreen tropical ecosystems. Several recent studies have demon- 
strated both the utility and challenges of remote sensing for studying 
tropical phenology (e.g., Huete et al., 2006; Myneni et al., 2007). 
The variability in integrated EVI for class 2 in this tile reflects this 





Fig. 8. Mean transition dates for different land cover classes in the C5 and C4 MLCD products. The mean values are computed for each land cover type present in each tile obtained 
from the C5 MODIS Land Cover Type product. 
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challenge and also suggests that MLCD results are probably less re- 
liable for tropical evergreen ecosystems. 

Finally, to illustrate the nature and magnitude of geographic 
patterns in interannual variability captured by the C5 MLCD product, 
Fig. 9 shows a map of anomalies in the timing of EVI increase and 
growing season length for 2002 relative to 2001-2006 averages 
(computed as 2002 minus the multi-year average), along with 
histograms for each. This figure suggests that the onset of green-up 
occurred later over much of North America relative to the 2001-2006 
average, especially at mid- to high latitudes and in the south-central 
United States. With the exception of the southeastern region, growing 
season anomalies follow the same general pattern and are positive 
(i.e., shorter growing season) throughout much of the continent. The 
climatic forcing behind this pattern is unclear, but it is likely that the 
widespread Northern Hemisphere drought that extended into 2002 
provides a partial explanation (Lotsch et al., 2005). 

3.4. Assessment using field data 

Validation of moderate resolution satellite products using ground 
measurements is challenging for several reasons including difficulties 
scaling plot level measurements to the resolution of the sensor field 
of view, geo-location uncertainties, limited temporal and spatial 
sampling of ground data, field instrument calibration, and sampling 
errors (Buermann et al., 2002; Herold et al., 2008; Huang et al., 2006; 
Tan et al., 2005a, b; Weiss et al., 2007). In addition to being relatively 
sparse, available field measurements related to phenology are 
normally collected for individual plants at scales well below the 


resolution of moderate spatial resolution sensor like MODIS. Further, 
the field of view of NBAR data generally includes a mosaic of species 
and landcover types (Baccini et al., 2008; Morisette et al., 2002). This 
presents significant challenges for product assessment and validation 
using measurements collected on the ground. 

In the specific case of the MLCD product, comparison between 
MODIS phenology and in-situ values is challenging for at least three 
reasons. First, the location of MODIS pixels and ground data may not 
match because of geo-location uncertainties (Tan et al., 2006). Second, 
the MODIS MLCD algorithm generates values for transition dates from 
NBAR EVI time series, which include uncertainties introduced by 
imperfect atmospheric correction, the presence of snow, cloud cover, 
and the 8-day temporal resolution of the EVI data. This issue is 
important because we are specifically interested in comparing events 
measured in time on the ground with estimated transition dates from 
the MLCD product, and best practices for achieving this goal are the 
subject of ongoing debate in the community (e.g., Schwartz & Hanes, 
2009; White et al, 2009). Third, in-situ measurements for a particular 
field site may not be representative of the associated 500-m MODIS 
pixel (Liang & Schwartz, 2009). Recognizing these limitations, here we 
compare MLCD results with available field data collected over multiple 
years from two long-term ecological reserve sites: (1) Harvard Forest 
in central Massachusetts and (2) Hubbard Brook in the White 
Mountain National Forest of New Hampshire. 

3.4.1. Comparison with Harvard Forest data 

Forest overstory and understory phenology data have been 
collected by Harvard Forest staff since 1989 using a consistent 
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Fig. 9. Anomalies in the timing of greenup onset and growing season length for 2002 relative to the 2001-2006 mean. 





S. Ganguly et al. / Remote Sensing of Environment 114 (2010) 1805-1816 


1813 


protocol and sampling strategy. Data are collected for two to five 
individuals of 33 native tree species located along a 2 km loop near the 
Harvard Forest headquarters (42°32'N, 72°11'W). Observations 
include measurements of percent canopy budburst and leaf develop- 
ment from April to June and percent leaf coloration and leaf drop from 
September to November at intervals of 3-7 days. Data from Harvard 
Forest were obtained from the Harvard Forest website (http://www. 
lternet.edu/hfr/). In this study, we examined data from 2001 to 2006. 

Table 3 presents a comparison of ground data collected at Harvard 
Forest with phenological transition dates estimated from MODIS for 
the onset timing of EVI increase (G 0 ), maximum (G M ), decrease (G s ), 
and minimum (G D ). Values for MODIS are based on a set of twenty- 
five C5 500-m MODIS pixels that intersect the ground data collection 
locations. For comparison, we include results for C4 (2001-2004). 
Corresponding indices based on ground measurements are also 
shown and defined in Table 3. Fig. 10 shows the time series of EVI 
data with mean values for phenological transition dates estimated by 
the MLCD algorithm at Harvard Forest in 2002. 

As we alluded to above, comparison of ground measurements of 
canopy phenology with MLCD transition dates is challenging. Here we 
compare continuous measures of canopy development (e.g., percent 
budburst) acquired at varying time intervals on the ground with 
MLCD transition dates. Specifically, we compare G 0 with the percent 
of the canopy for which budbreak had occurred on the same date; G M 
with the percent of the canopy that had achieved 75% and 95% of full 
leaf development; G s with the percent of the canopy with colored 
leaves and the percent of fallen leaves in the canopy; and G D with the 
average date of last green leaf in the canopy, and the percent of green 
leaves in the canopy at the time of G D . 

Values corresponding to the average budburst date (BBD) lag C5 
MODIS values of G 0 by about 1-17 days in each of the six years. In 
general G 0 values from MODIS correspond to the time when budburst 
has occurred for 5 to 33% of the forest canopy, and the range of BBD 
values in Table 3 includes the corresponding MODIS retrievals in all 
years except 2001 and 2002, when MODIS values for G 0 precedes the 
earliest non-zero value of BBD by 1 and 4 days, respectively. Overall, 
these results suggest that MODIS estimates of G 0 are sensitive to the 
early stages of leaf development at Harvard Forest. Note that MODIS 
results are also likely sensitive to understory phenology (not included 
in the data presented in Table 3), which tends to precede that of the 
forest canopy at Harvard Forest (Richardson & O’Keefe, 2009). 

The range of transition dates for the onset of EVI maximum (G M ) 
from MODIS corresponded to the timing when 87-100% of leaves had 
reached 75% of their final size (FL1), and when 72-97% of leaves had 
reached 95% their final size (FL2). These results suggest that the G M 



Fig. 10. Time series of mean EVI (±1 standard deviation) at Harvard forest for 2006 
with mean dates for the onset of EVI increase, maximum, decrease and minimum 
plotted as vertical lines. 

values correspond well to the time that maximum canopy develop- 
ment occurs. Note that G 0 and G M values in 2001 were probably 
affected by the MODIS instrument problems identified above. 

In the fall, average MODIS retrievals for the mean onset time of EVI 
decrease (G s ) range from DOY 224-241 over the six years considered. 
During this period the average percentage leaf coloring varied from 
0.3 to 1.8% and the average percentage of leaves that had fallen ranged 
from 0.1 to 0.8%. Average values for the onset time of EVI minimum 
(G d ) retrieved from MODIS differ from the average date of last green 
leaf coloring (AGD) by about 20 days. This result suggests that MLCD 
values for G D are biased late relative to the ground measurements. In 
general, however, comparison of spring and fall results with ground 
data suggest that MLCD algorithm is sensitive to seasonal transitions 
in EVI values and successfully captures phenological patterns at 
Harvard Forest (Table 3, Fig. 10). 

3.4.2. Comparison with Hubbard Brook field data 

The Hubbard Brook Experimental Forest is located within the 
White Mountain National Forest in central New Hampshire with 
elevations ranging from 222 to 1015 m above sea level. Researchers at 
Hubbard Brook have been monitoring the phenology of selected sugar 
maple, American beech and yellow birch trees at specific sites located 
in five watersheds since 1989 (Bailey, 2001). The phenological stage 


Table 3 

Comparisons among phenological transition dates retrieved from MLCD C4 and C5 data and in situ data collected at Harvard forest for the time period 2001-2006. For MODIS 
transition dates, the table presents the mean and standard deviation in DOY for the pixels included in the analysis. G 0 = onset date of EVI increase, G M = onset date of EVI maximum, 
G s = onset date of EVI decrease, and G D = onset date of EVI minimum; BBD = Date of first budburst in the field data; PC = percent canopy budburst on date of G 0 ; FL1 = Percent of 
canopy where leaves have reached 75% of final size; FL2 = Percent of canopy where leaves have reached 95% of final size; CL = percent of canopy with leaf coloring at time of G s ; 
FL = percent of canopy where leaves have fallen at time of G s ; AGD = average date of last green leaf in canopy; G L = percent of green leaves in the canopy at time of G D ; TL = percent 
of leaves remaining on the trees at time of mean G D . 


Year 

Go 




Gm 



G s 



G d 





MODIS 

BBD 



MODIS 

FL1 

FL2 

MODIS 

CL 

FL 

MODIS 

AGD 

GL 

TL 


X±(J 

x±cr 

Range 

PC 

x±cr 

x±cr 

x±cr 

X±(J 



X±(7 

X±C7 



2001 

108 ±5 (C5) 
119 ±2 (C4) 

125 ±9 

109-156 

8.25 

157 ± 7 (C5) 
158 ±3 (C4) 

87 ± 16 

72 ±23 

237 ±6 (C5) 

238 ±4 (C4) 

1 

0.16 

310±3 (C5) 
319 ±2 (C4) 

291 ±10 

0.91 

5.5 

2002 

102 ± 16 (C5) 
123 ± 1 (C4) 

117 ± 12 

106-153 

15.4 

173 ±7 (C5) 
167 ±2 (C4) 

98 ±2 

94 ±4 

241 ±11 (C5) 
233 ±3 (C4) 

0.3 

0.10 

315 ± 5 (C5) 
316±2 (C4) 

296 ±6 

0.00 

9 

2003 

125 ±4 (C5) 
129 ±2 (C4) 

126±9 

105-140 

33.0 

169 ±6 (C5) 
163 ±2 (C4) 

100 

90 ± 6 

233 ± 14 (C5) 

234 ±6 (C4) 

1.5 

0.16 

307 ±5 (C5) 
313 ±2 (C4) 

290 ±5 

0.02 

10 

2004 

111 ±4 (C5) 
119 ±2 (C4) 

123 ±6 

110-140 

5.0 

157 ±11 (C5) 
157±3 (C4) 

100 

87 ±7 

235 ±11 (C5) 

236 ±4 (C4) 

1.1 

0.10 

308 ±6 (C5) 
312 ±3 (C4) 

293 ±6 

0.18 

8 

2005 

114± 11 (C5) 

125 ± 10 

109-141 

16.5 

167 ±8 (C5) 

97 ±4 

88 ±8 

224 ± 12 (C5) 

0.9 

0.16 

315 ± 6 (C5) 

296 ±6 

0.00 

3 

2006 

115 ±3 (C5) 

123 ±7 

107-137 

15.0 

167 ±5 (C5) 

99 ± 1 

97 ±2 

230 ±4 (C5) 

1.8 

0.80 

303 ±5 (C5) 

287 ±9 

0.69 

10 
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of each tree is measured, through weekly visits in spring and autumn, 
on a canopy development index (CDI) that ranges from 0 to 4. Every 
plot record consists of observations on three trees each for sugar 
maple, yellow birch, and American beech. Values for the CDI depend 
on the season (spring or fall), and provide a measure of leaf 
emergence and canopy development in the spring, and leaf coloration 
and leaf drop in the fall. Specific details are available at http://www. 
hubbardbrook.org/data/dataset.php7id =51#. 

MLCD transition dates at Hubbard Brook and their associated 
standard deviations are listed in Table 4 along with canopy conditions 
quantified using the mean CDI at or close to MLCD transition dates. As 
for Harvard Forest, results from C4 are included for comparison. As we 
described above, ground measurement plots at Hubbard Brook span 
five watersheds. To allow comparisons with MLCD results, sixty-six 
500-m MODIS pixels covering the geographic extent of these water- 
sheds were averaged for this analysis. Also, because measurements at 
Hubbard Brook are collected weekly, it is difficult to compare these in- 
situ data with MODIS day of year transition dates. Hence, measured 
CDI values correspond to a fixed number of days (A) before or after 
the MODIS transition date. A A value of zero means that the reported 
canopy conditions coincide with the corresponding MODIS transition 
date. Negative values of A refer to the number of days that the in-situ 
measurements were collected before the MLCD transition date, and 
vice-versa. 

MLCD data are generally in good agreement with the ground 
measurements at Hubbard Brook. With the exception of G 0 in 2006, 
MLCD retrievals for G 0 , G M and G D broadly correspond to reasonable 
values of the CDI. However, values for G s indicate that canopy senes- 
cence begins in mid to late August, nearly 1 1 to 39 days before the CDI 
begins to decrease. Part of this discrepancy arises from the fact that in 
2001-2006 fall CDI data were not collected until early September. 
Further, documentation for the CDI indicates that a value of 4.0 in the fall 
allows for modest levels of leaf coloring. More generally, however, 
senescence tends to be quite gradual, and estimating a specific date 
associated with the onset of EVI decrease is difficult to detect from 
remote sensing. Thus, the results for G s probably reflect a combination of 
factors related to both data collection on the ground and lower overall 
quality in G s from MODIS relative to G 0 , G M and G D . Results for G s may 
also reflect the sensitivity of the EVI to seasonal variation in leaf 
properties such as water content (Cheng et al., 2006). 

As we noted above, MODIS G 0 values at Hubbard Brook in 2006 are 
unusually late with a mean value of 151, suggesting a problem with 
the algorithm or input data. Inspection of input data during this period 
reveals two related issues. First, NBAR snow flags indicate the pres- 
ence of snow at Hubbard Brook for the eight-day period starting on 
DOY 129. Second, NBAR blue band (MODIS band 3) data are missing 
during the same period, leading to missing EVI data during this period. 
The EVI time series for Hubbard Brook is shown in Fig. 11, where the 



Fig. 11. Time series showing NBAR- EVI data at Hubbard Brook before and after 
smoothing used to pre-process data prior to estimating transition dates. The data gap, 
filled by linear interpolation around DOY 129, is clearly evident. 

data gap is clearly evident. Thus, the combination of missing data and 
snow-flags in the input data causes the algorithm to spuriously 
estimate the onset of EVI increase well after DOY 129. Indeed, when 
the snow-flags from 2006 are replaced with values from 2005, the 
mean value for G 0 becomes DOY 121. This appears to be a relatively 
unusual situation, but requires correction in the next version of the 
product. 

4. Discussion and conclusions 

The goal of this paper was to present an overview and char- 
acterization of the Collection 5 MODIS land cover dynamics product 
through a series of inter-comparison exercises with the C4 product 
and available data from two field sites. In addition, the C5 product was 
evaluated to assess the nature and magnitude of variability captured 
by the product at interannual time scales. Results from a comparison 
of C4 and C5 data for six tiles spanning a range of climate and land 
cover types show that the C4 and C5 products are similar, with some 
notable differences. Most importantly, the C5 product has higher 
spatial resolution (500-m versus 1-km in C4) and fewer missing 
values outside of the tropics. C5 transition dates are generally con- 
sistent with C4, except in the tropics, where phenological dynamics 
are often more subtle and the frequency of missing data caused by 
clouds and aerosols is higher. At interannual time scales, observed 
year-to-year variability in C5 data is consistent with expected ranges 
associated with interannual climate variability. Larger differences at 


Table 4 

Mean and standard deviation in phenological transition dates (DOY) retrieved from C4 and C5 MODIS data at Hubbard Brook Experimental Forest for 2001-2006. A is the number of 
days before or after the MLCD value that in-situ data were collected, CDI is the canopy development index. Values in parentheses provide the range in CDI values on the day of year 
corresponding to the MLCD estimated transition date. 
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C D 




MODIS 
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CDI 

MODIS 

A 

CDI 

MODIS 

A 

CDI 

MODIS 

A 

CDI 


x±cr 



x±cr 



x±cr 



X ± O’ 



2001 

118 ±6 (C5) 
133 ± 6 (C4) 

10 

1.8 (1.4, 2.2) 

162 ±8 (C5) 
159±2 (C4) 

0 

4 (4,4) 

242 ±9 (C5) 
261 ±7 (C4) 

11 

3.5 (3.2, 3.7) 

299 ±8 (C5) 
306 ±2 (C4) 

-4 

0.3 (0,0.5) 

2002 

122 ± 10 (C5) 
138 ±5 (C4) 

4 

0.7 (0.3,1) 

169 ±7 (C5) 
168 ±5 (C4) 

-8 

4 (4,4) 

222 ±11 (C5) 
246 ±4 (C4) 

24 

4 (3.8,4) 

293 ±9 (C5) 
310±3 (C4) 

1 

1.4 (1,1.8) 

2003 

133 ±8 (C5) 
137 ±2 (C4) 

NA 

NA 

168 ±7 (C5) 
168 ±3 (C4) 

-1 

3.9 (3.8,4) 

226 ±9 (C5) 
241 ±4 (C4) 

32 

3.8 (3.6, 3.9) 

300 ±4 (C5) 
303 ±2 (C4) 

1 

0.3 (0.1, 0.4) 

2004 

116±6(C5) 
127±4 (C4) 

9 

1.1 (1,1.4) 

164 ±9 (C5) 
164 ±3 (C4) 

-11 

4 (4,4) 

226± 15 (C5) 
231 ±8 (C4) 

38 

3.7 (3.5, 3.9) 

299 ±6 (C5) 
304 ±3 (C4) 

0 

0.6 (0.3, 0.9) 

2005 

122 ± 11 (C5) 

0 

0.2 (0.1, 0.3) 

162 ±7 (C5) 

-5 

3.8 (3.7, 3.9) 

233 ±9 (C5) 

29 

4 (4,4) 

301 ±13 (C5) 

-4 

1.3 (0.9, 1.7) 

2006 

151 ± 6 (C5) 

-1 

3.6 (3.3, 3.8) 

173 ±13 (C5) 

-10 

4 (4,4) 

229 ±13 (C5) 

25 

3.8 (3.7, 3.9) 

289 ±5 (C5) 

0 

0.5 (0.3, 0.7) 
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local and regional scales may also reflect changing surface conditions 
in response to climate forcing, disturbance (fires and insect infesta- 
tions), and human management. Finally, MLCD transition dates 
from C5 were compared with field-measurements of forest canopy 
phenology at Harvard Forest and Hubbard Brook for 2001-2006. The 
results show that retrieved transition dates are generally realistic 
relative to field observations at each of these sites, but that MODIS- 
based estimates of phenological transition dates at the end of the 
growing season have larger uncertainty (and perhaps bias) relative to 
start-of-season metrics. 

Methods and data for assessment and validation of phenology 
metrics obtained from remote sensing is a critical issue requiring 
attention within the land surface phenology community. For this 
work, we relied on field data from two field sites in the northeastern 
United States. While these data provide measurements from a large 
number of trees and cover the period that MODIS has been collecting 
data, the limited geographic extent and ecological conditions in which 
these data were collected are clearly not sufficient to provide robust 
characterization of error in the MLCD product. Unfortunately, rela- 
tively few data sets are available that have been collected based on 
an explicit design in support of remote sensing product validation; 
this dearth of data is clearly limiting progress in this area (Liang & 
Schwartz, 2009). In particular, data sets collected across a wider range 
of ecosystems and that specifically measure phenology on the ground 
at spatial scales commensurate with moderate resolution remote 
sensing are urgently required. Intercomparison efforts such as those 
performed by White et al. (2009) provide useful guidance regarding 
differences and similarities among algorithm results. However, until 
data sets are available that have been collected using sample designs 
specifically implemented in support of remote sensing validation 
activities, it will be difficult to make meaningful conclusions regarding 
the strengths, weaknesses, and comparative accuracies of different 
algorithms and remote sensing products. Ongoing activities such as 
those coordinated by the National Phenology Network and other 
efforts using webcams (Richardson et al., 2009a, b) offer substantial 
promise for addressing this data gap, but current data is clearly not 
sufficient for global validation of land surface phenology products. 

The MODIS Land Cover Dynamics Product is one of a number of 
remote sensing-based products being used to generate regional to 
global scale maps of vegetation phenology. While the results 
presented in this and other work show that remote sensing can 
provide good quality results over large regions (e.g., temperate 
deciduous vegetation and agriculture), a number of important 
questions need to be resolved in order to both address weaknesses 
identified in this paper and to better meet the needs of scientists who 
wish to use these products. In addition to providing better 
characterization of the error and uncertainty associated with MLCD 
results, ongoing efforts are focused on developing improved methods 
for pre-processing input data (including screening for snow), and 
improved understanding of the nature and utility of retrieved 
phenological values in environments that present challenges for 
remote sensing including high latitude, arid, and tropical ecosystems. 
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